function Z = gauss2derror(c,data);

% See IniConditions2DFit() for 'c' structure
% c.amp - amplitude of 2D gaussian
% c.sigx - sigma_x
% c.sigy - sigma_y
% c.xoff - x-direction position offset of gaussian
% c.yoff - y-direction position offset of gaussian
% c.off - overall offset of image
% c.slopex - slope in the x-plane
% c.slopey - slope in the y-plane
% c.angle - angle of rotation of 2d gaussian

% Split the data matrix into x and y vectors
x = data(:, 1);
y = data(:, 2);

% Vector of errors
w = data(:, 3);

numbpoints = data(:,4);%matrix in which each entry is the number of points
sqrtnumbpoints = sqrt(numbpoints(1,1)); %sqrt of number of points

% FOR 2D GAUSSIAN WITH EQUALLY WEIGHTED DATA
Z =(c.off+c.slopex*(x-c.xoff)+slopey*(y-c.yoff)+c.amp*exp(-(((x-c.xoff).^2)/(2*c.sigx^2))-(((y-c.yoff).^2)/(2*c.sigy^2))));

% FOR ROTATED 2D GAUSSIAN WITHOUT SLOPES
%Z =offset+amp*exp(-(((x-xoff).*cosd(angle)+(y-yoff).*sind(angle))/(2*sigx)).^2-((-(x-xoff).*sind(angle)+(y-yoff)*cosd(angle))/(2*sigy)).^2);

% FOR 2D GAUSSIAN WITH WEIGHTED DATA
%Z =(offset+slopex*(x-xoff)+slopey*(y-yoff)+amp*exp(-(((x-xoff).^2)/(2*sigx^2))-(((y-yoff).^2)/(2*sigy^2))))./(w*sqrtnumbpoints); % We're trying to fit z = f(x, y) so compute f(x, y)
